Use R package UScensus2010county to complete the following tasks: (20 pt.)
library(UScensus2010county)
Plot a map of New York counties using the plot function.
data(new_york.county10)
shp <- new_york.county10
df <- shp@data
df
plot(shp)
Plot a map of New York counties using the qtm function.
install.packages("tmap")
Error in install.packages : Updating loaded packages
library(tmap)
qtm(shp)
How many counties in New York State?
df
print("There are 62 counties in New York State")
[1] "There are 62 counties in New York State"
What’s the 3-digit fips code of Broome County?
print("36007")
[1] "36007"
Compute descriptive statistics of the population column (P0010001), including total, minimum, maximum, mean, median, and skewness.
nypop <- df$P0010001
t <- sum(nypop)
t
[1] 19378102
a <- summary(nypop)
a
Min. 1st Qu. Median Mean 3rd Qu. Max.
4836 51244 91301 312550 231060 2504700
library(moments)
b <- skewness(nypop)
b
[1] 2.579801
Create a histogram and a boxplot of the population.
a <- hist(nypop)
a
$breaks
[1] 0 500000 1000000 1500000 2000000 2500000 3000000
$counts
[1] 53 3 3 1 1 1
$density
[1] 1.709677e-06 9.677419e-08 9.677419e-08 3.225806e-08 3.225806e-08 3.225806e-08
$mids
[1] 250000 750000 1250000 1750000 2250000 2750000
$xname
[1] "nypop"
$equidist
[1] TRUE
attr(,"class")
[1] "histogram"
b <- boxplot(nypop)
b
$stats
[,1]
[1,] 4836
[2,] 51125
[3,] 91301
[4,] 234878
[5,] 468730
$n
[1] 62
$conf
[,1]
[1,] 54429.09
[2,] 128172.91
$out
[1] 744344 1585873 2504700 1339532 949113 2230722 919040 1493350 1385108
$group
[1] 1 1 1 1 1 1 1 1 1
$names
[1] "1"
Use R package UScensus2010tract to complete the following tasks: (20 pt.)
library(UScensus2010)
if(!require(UScensus2010tract)) install.tract("osx")
library(UScensus2010tract)
Plot a map of New York census tracts using the plot function.
data(new_york.tract10)
shp <- new_york.tract10
df <- shp@data
df
Restarting R session...
Compute the total population based on census tracts.
x
[1] 19378102
Select all census tracts in Broome County and plot the map.
BroomeCounty <- county(fips="007",name="Broome",state="NY",level="tract")
plot(BroomeCounty)
What’s the total population of Broome County?
Broomepop <- sum(BroomeCounty$P0010001)
Broomepop
[1] 200600
Create a histogram and a boxplot of population based on census tracts of Broome County.
hist(BroomeCounty$P0010001)
boxplot(BroomeCounty$P0010001)
Select the first five columns of the shapefile of Broome County census tract; add a new population ratio column (= census tract population / county population); save the new shapefile to the hard drive.
pr <- data.frame(BroomeCounty[,1:5])
pr$ratio <- BroomeCounty$P0010001/200600
pr
Use R packages raster and ncdf4 to complete the following tasks: (20 pt.)
library(raster)
library(ncdf4)
Load the multi-band raster dataset “NDVI.nc” into R.
setwd("data/gimms_ndvi/")
Error in setwd("data/gimms_ndvi/") : cannot change working directory
Get the basic information about the dataset, including the number of rows, columns, cells, and bands; spatial resolution, extent, bounding box, and projection.
ndvi.rb <- brick("NDVI.nc")
ndvi.rb
class : RasterBrick
dimensions : 360, 720, 259200, 36 (nrow, ncol, ncell, nlayers)
resolution : 0.5, 0.5 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
data source : C:\STATS\Lab_10\NDVI.nc
names : X2000.01.01, X2000.02.01, X2000.03.01, X2000.04.01, X2000.05.01, X2000.06.01, X2000.07.01, X2000.08.01, X2000.09.01, X2000.10.01, X2000.11.01, X2000.12.01, X2001.01.01, X2001.02.01, X2001.03.01, ...
Date : 2000-01-01, 2002-12-01 (min, max)
varname : NDVI
Aggregate all bands to generate a mean NDVI raster and a maximum NDVI raster; save the two new raster datasets to the hard drive.
mean(ndvi.rb)
class : RasterLayer
dimensions : 360, 720, 259200 (nrow, ncol, ncell)
resolution : 0.5, 0.5 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
data source : in memory
names : layer
values : 0, 0.8537528 (min, max)
writeRaster(mean(ndvi.rb),filename = "mean_ndvi.tif", overwrite=TRUE)
max(ndvi.rb)
class : RasterLayer
dimensions : 360, 720, 259200 (nrow, ncol, ncell)
resolution : 0.5, 0.5 (x, y)
extent : -180, 180, -90, 90 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
data source : in memory
names : layer
values : 0, 0.9681 (min, max)
writeRaster(max(ndvi.rb),filename = "max_ndvi.tif", overwrite=TRUE)
Plot the maps of monthly NDVI of the year 2001
Create histograms of monthly NDVI of the year 2001.
hist(ndvi2001)
size changed to n because it cannot be larger than n when replace is FALSE39% of the raster cells were used (of which 42% were NA). 58093 values used.
Plot the NDVI map of July 2000; click any location with data on the map and retrieve the NDVI time series for all years; plot the NDVI time series of the selected location.
plot(ndvi.rb,7)
values <- ndvi.rb[150, 50]
values <- click(ndvi.rb, n=1, xy=TRUE)
values <- click(ndvi.rb, n=1, xy=FALSE)